RVAideMemoire
パッケージにあるfisher.multcomp
関数で,p.adjust
関数にあるp値調整法でFisher検定の多重比較ができる.
library(RVAideMemoire)
<- matrix(c(17,23,12,24,20,10),ncol=2,
d dimnames=list(c("cont","trt1","trt2"),
c("alive","dead")),
byrow=TRUE)
d
## alive dead
## cont 17 23
## trt1 12 24
## trt2 20 10
fisher.test
関数データ例で,cont
,trt1
,trt2
の3群の生存率の比較をFisher検定で行う.
<- fisher.test(d[c(1,2),])
res12 <- fisher.test(d[c(1,3),])
res13 <- fisher.test(d[c(2,3),])
res23 c(`cont vs trt1`=res12$p.value,
`cont vs trt2`=res13$p.value,
`trt1 vs trt2`=res23$p.value)
## cont vs trt1 cont vs trt2 trt1 vs trt2
## 0.48195094 0.05567911 0.01279792
RVAideMemoire
パッケージのfisher.multcomp
関数で, p.method="none"
で調整なしのFisher検定のp値.
fisher.multcomp(d, p.method = "none")
##
## Pairwise comparisons using Fisher's exact test for count data
##
## data: d
##
## cont trt1
## trt1 0.48195 -
## trt2 0.05568 0.0128
##
## P value adjustment method: none
\(m\)個の帰無仮説\(H_1,H_2,\ldots,H_m\)の検定をする状況を考える. \(m_0\)個の帰無仮説には真に正しい仮設とする. 下の表は,\(V\)個が第1種の過誤(真である帰無仮説を誤って棄却)が生じた検定数であり, \(R\)は棄却した検定数である. \(R\)は観測可能な確率変数に対し,\(S\), \(T\), \(U\), \(V\)は観測できない確率変数である. \(m\)と\(m_0\)は定数であり,\(m_0\)は未知である.
\[\begin{array}{cccc} \hline {\rm Hypothesus} & {\rm Not\ Rejected} & {\rm Rejected} & {\rm Total} \\ \hline {\rm True} & U & V & m_0 \\ {\rm False} & T & S & m-m_0 \\ \hline {\rm Total} & m-R & R & m \\ \hline \end{array}\]個々の帰無仮説\(H_i\)に対する検定における第1種の過誤は有意水準\(\alpha\)で制御されていが, \(m\)個の帰無仮説の検定の全体での第1種の過誤は,FWER(Family Wise Error Rate)で考える.
\[ {\rm FWER}=\Pr(V>0) \]
次に,棄却した検定数\(R\)に含まれる偽陽性(False Positive)の検定数として\(V\)を扱う方法を考える. 間違って棄却された帰無仮説の割合を\(Q=V/R\)とする.ただし,\(R=0\)のときは\(Q=0\)とする. この\(Q\)も確率変数である.その期待値をFDR(False Discovery Rate)と呼ぶ. \[ {\rm FDR}=E(Q)=E\left(\frac{V}{R} \right) \]
FWERとFDRには以下の関係がある.
\(m\)個の検定の\(p\)値を値が小さい順に並べ替えた\(p_1<p_2< \cdots <p_m\)に対して, ボンフェローニ法では,調整\(p\)値を\(q_i=\min(1,mp_i)\) \((i=1,2,\ldots,m)\)とする. 上の例だと,\(p_1=0.0128\),\(p_2=0.05568\),\(p_3=0.48195\)より,
\[q_1=3\times 0.0128=0.0384, \quad p_2=3\times 0.05568=0.1670, \quad p_3=3\times0.48195 \to 1\]
fisher.multcomp
関数でp.method="bonferroni"
とする.
fisher.multcomp(d, p.method = "bonferroni")
##
## Pairwise comparisons using Fisher's exact test for count data
##
## data: d
##
## cont trt1
## trt1 1.000 -
## trt2 0.167 0.03839
##
## P value adjustment method: bonferroni
ホルム法では次のように\(p\)値を調整し,\(q_i=(m-k+1)p_{k}<\alpha\)となる最大の\(k\)を見つけ,第1順位から第\(k\)順位までの帰無仮説を棄却し,第\((k+1)\)順位以降の帰無仮説の判定を保留する.
\[ q_1=mp_1,\ q_2=(m-1)p_2,\ q_3=(m-2)p_2,\ \cdots \ , q_{m-1}=2p_{m-1},\ q_m=p_1,\ \]
上の例だと,\(p_1=0.0128\),\(p_2=0.05568\),\(p_3=0.48195\)より,
\[q_1=3\times 0.0128=0.0384, \quad p_2=2\times 0.05568=0.11136, \quad p_3=1\times 0.48195=0.48195\]
fisher.multcomp
関数でp.method="holm"
とする.
fisher.multcomp(d, p.method = "holm")
##
## Pairwise comparisons using Fisher's exact test for count data
##
## data: d
##
## cont trt1
## trt1 0.4820 -
## trt2 0.1114 0.03839
##
## P value adjustment method: holm
Benjamini-Hochberg法では次のように\(p\)値を調整し,\(q_i=(m/i)p_i < \alpha\)となる最大の\(k\)を見つけ, \(H_1,H_2,\ldots,H_k\)の仮設を棄却する. \(q_1<\alpha\)を満たさなかった場合は,どの帰無仮説も棄却しない.
\[ q_1=\frac{m}{1}p_1,\ q_2=\frac{m}{2}p_2,\ q_3=\frac{m}{3}p_2,\ \cdots \ , q_{m-1}=\frac{m}{m-1}p_{m-1},\ q_m=\frac{m}{m}p_m,\ \]
fisher.multcomp(d, p.method = "BH")
##
## Pairwise comparisons using Fisher's exact test for count data
##
## data: d
##
## cont trt1
## trt1 0.48195 -
## trt2 0.08352 0.03839
##
## P value adjustment method: BH